import numpy as np
import save_contour
save_contour = reload( save_contour )
import __main__
As11t,As11c,As22t,As22c,As12,As13 = __main__.cc.allowables 

def get_allow( s11, s22 ):
    global As11c, As11t, As22c, As22t, As12
    if   s11 <  0 and s22 <  0:
        return As11c, As22c, As12
    elif s11 <  0 and s22 >= 0:
        return As11c, As22t, As12
    elif s11 >= 0 and s22 <  0:
        return As11t, As22c, As12
    elif s11 >= 0 and s22 >= 0:
        return As11t, As22t, As12

def nap( cc, pos ):
    num   = int(pos/3) +1
    angle = cc.stack[ num-1 ]
    tmp   = ['bot','mid','top']
    pos   = tmp[ int(pos) % 3 ]
    return num, angle, pos

def envelope_values( values_dict, positions_dict, minmax ):
    ans     = dict([ [k,None] for k in values_dict.keys() ])
    ans_pos = dict([ [k,None] for k in values_dict.keys() ])
    for key in values_dict.keys():
        values = values_dict[ key ].values
        positions = positions_dict[ key ].values
        if   minmax == 'min':
            minvalue     =  1.e6
            minvalue_pos = -1
            for i in range(len(values)):
                v   = values[ i ]
                if minvalue > v.data:
                    minvalue     = v.data
                    pos = positions[ i ]
                    minvalue_pos = pos.data
            ans[key]     = minvalue
            ans_pos[key] = minvalue_pos
        elif minmax == 'max':
            maxvalue     = -1.e6
            maxvalue_pos = -1
            for i in range(len(values)):
                v   = values[ i ]
                if maxvalue < v.data:
                    maxvalue     = v.data
                    pos = positions[ i ]
                    maxvalue_pos = pos.data
            ans[key]     = maxvalue
            ans_pos[key] = maxvalue_pos
    return ans, ans_pos

def calc_frame( self, frame, check_print_report = True ):
    global As11t, As11c, As22t, As22c
    f = frame
    total = len(self.zdisp) - 1
    print 'processing frame % 4d / % 4d, zdisp %1.3f, zload, %1.2f' \
          % (f.frameId, total, self.zdisp[ f.frameId ], self.zload[ f.frameId ])
    hashin         = {'HSNFCCRT':None, 'HSNFTCRT':None, 'HSNMCCRT':None, 'HSNMTCRT':None}
    hashin_max     = {'HSNFCCRT':None, 'HSNFTCRT':None, 'HSNMCCRT':None, 'HSNMTCRT':None}
    hashin_max_ms  = {'HSNFCCRT':None, 'HSNFTCRT':None, 'HSNMCCRT':None, 'HSNMTCRT':None}
    hashin_max_pos = {'HSNFCCRT':None, 'HSNFTCRT':None, 'HSNMCCRT':None, 'HSNMTCRT':None}
    hashin_pts     = {'HSNFCCRT':None, 'HSNFTCRT':None, 'HSNMCCRT':None, 'HSNMTCRT':None}
    hashin_fields  = {'HSNFCCRT':[]  , 'HSNFTCRT':[]  , 'HSNMCCRT':[]  , 'HSNMTCRT':[]  }
    stress         = {'S11':None, 'S22':None, 'S12':None }
    stress_min     = {'S11':None, 'S22':None, 'S12':None }
    stress_min_ms  = {'S11':None, 'S22':None, 'S12':None }
    stress_min_pos = {'S11':None, 'S22':None, 'S12':None }
    stress_max     = {'S11':None, 'S22':None, 'S12':None }
    stress_max_ms  = {'S11':None, 'S22':None, 'S12':None }
    stress_max_pos = {'S11':None, 'S22':None, 'S12':None }
    stress_pts     = {'S11':None, 'S22':None, 'S12':None }
    stress_fields  = {'S11':[]  , 'S22':[]  , 'S12':[]   }
    for key in hashin.keys():
        hashin[key]     = f.fieldOutputs[key]
        hashin_pts[key] = f.fieldOutputs[key].locations[0].sectionPoints
    for key in stress.keys():
        stress[key]     = f.fieldOutputs['S'].getScalarField(componentLabel=key)
        stress_pts[key] = f.fieldOutputs['S'].locations[0].sectionPoints
    import abaqus
    for key in hashin.keys():
        for pt in hashin_pts[key]:
            hashin_fields[key].append( hashin[key].getSubset( sectionPoint=pt ))
        hashin_max[key], hashin_max_pos[key] =\
            abaqus.maxEnvelope( hashin_fields[key] ) 
    for key in stress.keys():
        for pt in stress_pts[key]:
            stress_fields[key].append( stress[key].getSubset( sectionPoint=pt ))
        stress_min[key], stress_min_pos[key] = \
            abaqus.minEnvelope( stress_fields[key] )
        stress_max[key], stress_max_pos[key] = \
            abaqus.maxEnvelope( stress_fields[key] )
    stress_min_num, stress_min_pos_num =\
        envelope_values( stress_min, stress_min_pos, 'min')
    stress_max_num, stress_max_pos_num =\
        envelope_values( stress_max, stress_max_pos, 'max')
    hashin_max_num, hashin_max_pos_num =\
        envelope_values( hashin_max, hashin_max_pos, 'max')
    As11, As22, As12 = get_allow( stress_min_num['S11'], stress_min_num['S22'] )
    stress_min_ms['S11'] = abs(As11 / stress_min_num['S11']) - 1.
    stress_min_ms['S22'] = abs(As22 / stress_min_num['S22']) - 1.
    stress_min_ms['S12'] = abs(As12 / stress_min_num['S12']) - 1.
    As11, As22, As12 = get_allow( stress_max_num['S11'], stress_max_num['S22'] )
    stress_max_ms['S11'] = abs(As11 / stress_max_num['S11']) - 1.
    stress_max_ms['S22'] = abs(As22 / stress_max_num['S22']) - 1.
    stress_max_ms['S12'] = abs(As12 / stress_max_num['S12']) - 1.
    for key in hashin_max_num.keys():
        hashin_max_ms[key] = 1./np.sqrt( hashin_max_num[key] ) - 1.  
    self.stress_min_num[frame.frameId]     = stress_min_num
    self.stress_min_ms[frame.frameId]      = stress_min_ms
    self.stress_min_pos_num[frame.frameId] = stress_min_pos_num
    self.stress_max_num[frame.frameId]     = stress_max_num
    self.stress_max_ms[frame.frameId]      = stress_max_ms
    self.stress_max_pos_num[frame.frameId] = stress_max_pos_num
    self.hashin_max_num[frame.frameId]     = hashin_max_num
    self.hashin_max_ms[frame.frameId]      = hashin_max_ms
    self.hashin_max_pos_num[frame.frameId] = hashin_max_pos_num
    if check_print_report:
        print '\t'.join([ '','Hashin,FC','Hashin,FT','Hashin,MC','Hashin,MT',
                       'S11min','S11max','S22min','S22max','S12min','S12max'])
        hsn  = self.hashin_max_num[ frame.frameId ]
        smin = self.stress_min_num[ frame.frameId ]
        smax = self.stress_max_num[ frame.frameId ]

        print '\t'.join([ str(i) for i in [
          'stress',
          hsn['HSNFCCRT'], hsn['HSNFTCRT'], hsn['HSNMCCRT'], hsn['HSNMTCRT'],
          smin['S11'],     smax['S11'],     smin['S22'],     smax['S22'],
          smin['S12'],     smax['S12']
                                          ]])
        hsnpos  = self.hashin_max_pos_num[ frame.frameId ]
        sminpos = self.stress_min_pos_num[ frame.frameId ]
        smaxpos = self.stress_max_pos_num[ frame.frameId ]
        indexes = [0,1,2]
        names   = ['ply num', 'ply angle', 'ply pos']
        for i in indexes:
            print '\t'.join([ str(i) for i in [
              names[i],
              nap(self,hsnpos['HSNFCCRT'] )[i],nap(self,hsnpos['HSNFTCRT'])[i],
              nap(self,hsnpos['HSNMCCRT'] )[i],nap(self,hsnpos['HSNMTCRT'])[i],
              nap(self,    sminpos['S11'] )[i],nap(self,    smaxpos['S11'])[i],
              nap(self,    sminpos['S22'] )[i],nap(self,    smaxpos['S22'])[i],
              nap(self,    sminpos['S12'] )[i],nap(self,    smaxpos['S12'])[i],
                                              ]])
        print '\t'.join([ str(i) for i in [
            'allowables', 1.,1.,1.,1.,As11c,As11t,As22c,As22t,As12,As12
            ]])
        hsnms  = self.hashin_max_ms[ frame.frameId ]
        sminms = self.stress_min_ms[ frame.frameId ]
        smaxms = self.stress_max_ms[ frame.frameId ]
        print '\t'.join([ str(i) for i in [
            'margin of safety',
            hsnms['HSNFCCRT'], hsnms['HSNFTCRT'],
            hsnms['HSNMCCRT'], hsnms['HSNMTCRT'],
            sminms['S11']    , smaxms['S11']    ,
            sminms['S22']    , smaxms['S22']    ,
            sminms['S12']    , smaxms['S12']
                                          ]])

def calc_frames( self, frames=None, MSlimits=[0.0,0.2,0.5], frame_indexes=[] ):
    if frames == None:
        frames = self.attach_results().steps[ self.step2Name ].frames
    if frame_indexes == []:
        iframes = frames
    else:
        iframes = [ frames[i] for i in frame_indexes ]
    print_ms_tsai_hill = [True for i in MSlimits]
    print_ms_tsai_wu   = [True for i in MSlimits]
    i_mslimit = 0
    # finding frame first buckling and first minimum
    fb_load = 0.
    for i in range(len(iframes)):
        zload = self.zload[i]
        if abs(zload) < abs(fb_load):
            break
        fb_load = zload
    fb_frame = iframes[i-1]
    fb_load  = self.zload[ fb_frame.frameId ]
    for i in range(len(iframes)):
        if i <= fb_frame.frameId:
            continue
        zload = self.zload[i]
        if abs(zload) > abs(fb_load):
            break
        fb_load = zload
    fm_frame = iframes[i-1]                
    #
    print '\nreport at first buckling load'
    calc_frame( self=self, frame=fb_frame, check_print_report=True )
    print '\nreport at first minimum'
    # calculating minimum required stiffness
    delta_u = self.zdisp[ fm_frame.frameId ] - self.zdisp[ fb_frame.frameId ]
    delta_f = self.zload[ fm_frame.frameId ] - self.zload[ fb_frame.frameId ]
    k_min   = abs( delta_f / delta_u ) 
    print 'elastic recovery, minimum stiffness %1.3f' % k_min
    calc_frame( self=self, frame=fm_frame, check_print_report=True )
    #
    self.stress_min_num     = {}
    self.stress_min_ms      = {}
    self.stress_min_pos_num = {}
    self.stress_max_num     = {}
    self.stress_max_ms      = {}
    self.stress_max_pos_num = {}
    self.hashin_max_num     = {}
    self.hashin_max_ms      = {}
    self.hashin_max_pos_num = {}
    for f in iframes:
        calc_frame( self=self, frame=f, check_print_report=False )
    # saving contours
    #save_contour.save_contour( self=self, frame_num=fb_frame.frameId,
    #                           sufix='first_buckling')
    #save_contour.save_contour( self=self, frame_num=fm_frame.frameId,
    #                           sufix='first_minimum' )
#    for i in range( len(report_frames) ):
#        frame = report_frames[ i ]
#        sufix = sufixes[ i ]
#        save_contour.save_contour( self=self, frame_num=frame.frameId, sufix=sufix )
#        print '\nreport at MS %1.2f', MSlimits[ i ] 
#        # calculating minimum required stiffness
#        delta_u = self.zdisp[ frame.frameId ] - self.zdisp[ fb_frame.frameId ]
#        delta_f = self.zload[ frame.frameId ] - self.zload[ fb_frame.frameId ]
#        k_min   = abs( delta_f / delta_u ) 
#        print 'elastic recovery, minimum stiffness %1.3f' % k_min
#        calc_frame( self=self, frame=frame, check_print_report=True )
